Elongated Fermi superfluid: absence of critical imbalance enhancement at equilibrium 
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We show that the maximum population imbalance ratio Pec for a two-component Fermi gas near the unitarity 
limit to condense does not increase with the trap aspect ratio A, by two methods of 1) solving the Bogoliubov- 
de Gennes equations with coupling-constant renormalization, and 2) studying the pairing susceptibility by the 
real-space self-consistent T-matrix approximation. The deviation of the cloud shape from what is expected 
from the trap shape increases but stays minor with increasing A up to 50. This finding indicates that despite the 
apparent discrepancy between the MIT and Rice experiments over the value of Pec and the validity of local 
density approximation, the equilibrium state of the system for the aspect ratio in the Rice experiment should be 
consistent with that of MIT. 
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Gaseous Fermi superfluids are endowed with new de- 
grees of controllability over population difference and trap 
anisotropy. Imbalanced superfluidity of ^Li has been ob- 
served by the Rice HI] and MIT groups, but their results 
have shown marked differences over the validity of local den- 
sity approximation (LDA) and the Chandrasekhar-Clogston 
(CC) limit - the upper bound Pqc of imbalance parameter 
P = {N^ — Ni)/N beyond which superfluidity breaks down 
ijstl, where N-^- and Ni are the numbers of majority and mi- 
nority atoms, and N = + Ni is the total atom number. 
In the MIT experiment the profiles of both majority and mi- 
nority clouds obey LDA, while in the Rice experiment with a 
very elongated trap and fewer atoms, LDA apparently breaks 
down. The CC limit was observed at MIT but not at Rice. 
A phenomenological surface tension [4,^ of the condensate 
was shown to reproduce the deformation observed by the Rice 
group, but how to reconcile the apparently contradicting ex- 
perimental differences without free parameters remains elu- 
sive . More recently the non-equilibrium state during the 
evaporative cooling process ['s'l was discussed to explain the 
Rice results. In this Letter we demonstrate, for the equilib- 
rium state of the system at low and finite temperatures, that 1) 
the CC limit does not increase with increasing the trap aspect 
ratio A, and that 2) while the density-difference distribution 
does deform from what is expected from the trap shape, the 
deformation is not as significant as in the Rice experiment for 
the number of atoms as small as 3 x 10*. 

We consider a system of atoms with mass m confined in 
an axisymmetric harmonic potential V{r) = m(cj^(a;'^ + 
y^) +w^z^)/2 with axial frequency Uz and radial one uj±, and 
analyze superfluidity of this system using the Bogoliubov-de 
Gennes (BdG) equations [9-16]. Sensarma ef a/. [15] studied 
the shape of the atom cloud by changing N and P{< 0.4), 
and ai-gued that (iV/A)^/^ > 1 (A = uj±/uJz) should be the 
condition for the vahdity of LDA. For (7V/A)^/^ 10, the 
cloud shape obtained in Ref. [15] looks quite similar to that 
of the equipotential surface. However, the Rice experiment 



shows the breakdown of LDA for almost the same value of 
(AT/A) 1/3. While our numerical results also show deformation 
similar to that found in Ref. [15J for TV ~ 10^ at A 4, the 
density profiles are different presumably because we incorpo- 
rate the effect of the chemical potential difference as well as 
the interaction between atoms in the normal state. Such defor- 
mation almost disappears for ^ 3 x 10**. 

The BdG equations for unequal chemical potentials 
(Mt'Mi) are given by 



A* -Hi - W-f 



(1) 



where = - V^/(2m) + V{r) - fi^ (ct =t, i) is the one- 
body Hamiltonian, and Wa{r ) is the Hartree-Fock mean-field 
energy gna{r) with the coupling constant g given in terms of 
s-wave scattering length Os as g ~ Airh^as/m. In the follow- 
ing we take to = ft = fee = 1, set ZU = y/U^Iol — ijJ_l/ v^A, 
and choose y^h/{muJ) = 1 as the unit of length. The self- 
consistent conditions give the density distributions ncr{r) and 
the s-wave singlet pair amplitude A(r) as 

Mr) = E fiMr)\', %(r) = ^ (1 - /,) Mr)\' , 

A{r)=g,s{r)J2f,Mr)v;{r), (2) 

Q 

where fq = (e^'^i + 1)^^ is the Fermi distribution function 
with 13 = {kBT)~^. To cope with the ultraviolet divergence in 
A(r), we follow Bulgac and Yu [ 17] and treat the contribution 
from states above an energy cutoff Ec within LDA. In Ref. 
HtIi . where /i = /i^ = is assumed, the single-particle 
Green's function G° with Hq = - V^/(2to) + F - /i is used 
to remove the divergence. The regular part G^'"^"^ of Gj^ is 
obtained by employing the Thomas-Fermi approximation for 
the states above Ec, so that the effective coupling constant 
is given in terms of fcc(r) = •\/2 {Ec — V{r)) and fcp(r) = 
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FIG. 1. (Color online) Top panel shows the P dependence of the ra- 
tio R of the axial to radial cloud widths at which the atomic density 
equals 1/20 of its peak value. Lower panels show the density differ- 
ence between majority and minority atoms plotted for A = 5 (left) 
and 50 (right) with fcpas = —1.3, N = 3 x 10*, and varying imbal- 
ance parameter P. The trap axis lies in the horizontal direction. At 
A = 5, (Mt>M4) = (38.83,28.80) for P = 0.1 and (53.08,4.79) 
for P = 0.8; at A = 50, {fi^, fii) = (38.63, 28.78) for P = 0.1 and 
(53.13,4.66) forP = 0.8. 
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Grasso and Urban ifisll replaced fcp(r) with kpir) 
^y2{^-V{r) ~W{r)), where W{r) = W^^i for /i^ = /z^, 
so that the convergence is achieved for much smaller values 
of Ec- We adopt this method except that we replace Gj^ by 

^G^'^'^'s + G°'^'^°s^ /2 to maintain a given chemical potential 
difference. Consequently, Eq. (O is replaced by 



(4) 
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where kp^ir) = y/2 {fi^ ~ V{r) - WV(r)), with W^r) = 
gn^ir). While BdG theory was originally proposed to de- 
scribe the weak-coupling BCS limit, it was demonstrated to 



describe the BEC limit fl?], and the BCS-BEC crossover re- 
gion was also studied by this theory [20]. We therefore ex- 
pect that this theory is applicable, at least qualitatively, for the 
strongly interacting region with population imbalance, pro- 
vided that an appropriate coupling-constant renormalization 
is employed. 

At the unitarity limit (/cpas)^^ — > 0, the normal state in- 
teraction does not diverge, and the binding energy of a single 
I atom to the Fermi sea of ^ atoms with the Fermi energy 
Spt is -(3/5)yl£;Ft with A = 0.97(2) [21]. This coitc- 
sponds to the mean-field energy of — (9yl7r/20) {kY^{r))^^ x 
47rn^(r)n4.(r), where k-pa = (eTr^n^)^/^. On the BCS 
side of the unitarity limit, the normal state interaction should 
be weaker than at the unitarity limit so |fcF-f-(r)as| < 
|fcF(0)as| < 9A7r/20 = 1.37. Moreover, we can show that 
the BdG equations ([T|) do not have a stable self-consistent so- 
lution for Ifepflsl > 37r/4 = 2.36 for the homogeneous case 
without chemical potential difference. 

We take k^T — 0.05SaJ and use the Steffensen iteration to 
solve Eqs. ([T]l and to self-consistently determine n^{r), 
ni,{r) and A(r) for a given set of (/i-|-, /ij,). The number of 
atoms in the cr{=^,\) state is defined as = J (Prn„. 

Figure [T] shows the main results of this Letter. For both 
A = 5 and A = 50, the ratio of the axial to radial cloud widths 
remains close to A for both minority and majority atoms, and 
the dip of the density difference rapidly dwindles with increas- 
ing P, vanishing for P > 0.75. Thus the CC limit is not 
enhanced as A is increased. For A = 50, the density differ- 
ence shows some deformation for small P, but it disappears 
for P > 0.60. 

Figure |2] shows typical distributions of n^^i, their differ- 
ence, and A for P ~ 0.40. We rescale the calculated distri- 
bution as r ^ A?' so that the equipotential surface becomes a 
circle. For A = 5, the shape of the minority component and 
the density difference closely follow the equipotential surface, 
as shown in the left column of Fig. [T] The pair amplitude 
shows sign changes, which are absent in LDA but shows up 
in the BdG simulation as discussed in Ref. Ii 1 3il for a spherical 
system. 

For P larger than 0.7, the pair amplitude almost vanishes, 
and the density difference peaks at r — 0. We therefore con- 
clude that LDA is essentially valid at A = 5 as observed by the 
MIT group. For A = 50, while the density difference shows 
some deviation from the trap shape, implying the breakdown 
of LDA, the degree of breakdown is rather small. This can be 
seen from almost spherical density distributions of both the 
majority and minority components in Fig. |2] (Note that in 
Fig.|2]the vertical axis is scaled by a factor of A.) The region 
with non-vanishing pairing amplitude is also rather similar to 
that of the minority component, reflecting the fact that pairing 
occurs effectively in the strongly interacting regime. 

With the same number of atoms, we have thus confirmed 
that LDA is less invalid at A = 5 than at A = 50. The 
breakdown of LDA is a finite-size effect, and it is enhanced 
for larger A. Figure |3] shows the atom-number dependence 
of Pec for a spherical trap. We find that with increasing N, 
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FIG. 2. (Color online) Majority and minority density distributions n^^i{z,r), their difference n^{z,r) — ni{z,r) and pair amplitude A{z,r) 
plotted for (a) A = 5 and (b) A = 50 with k^tas = -1.3 and iV = 3 x 10** at P = 0.40. (^t, Mi) = (45.90, 19.20) in (a) and (45.57, 19.42) 
in (b). The density distributions and A are displayed in color-coded gauges shown on the upper right and lower right, respectively. In each 
inset, the cross sections at r = (solid curve) and z — (dashed curve) are plotted against z and Xr, respectively. For the pair amplitude, the 
regions close to the horizontal axis are enlarged in the smaller insets. 
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FIG. 3. (Color online) CC limit Pec plotted against the total number 
of trapped atoms A'^ in a spherical harmonic potential. Here, Pec is 
identified as the value of P at which the extrapolated plot of A(0) 
crosses zero. 



FIG. 4. (Color online) Largest eigenvalue AMax of xsc obtained in 
the RSTA method is plotted against T/ [ji^/k^) for = lOfiZZJ = 
lQh\/u)zU}'^ for aspect ratios A = 1 (open symbols) and 12 (filled) 
and imbalance parameters P = 0.1, 0.2, 0.3, 0.4. 



Pec approaches a constant value for each jfep-f^asl, which is, 
for l/cF-t-flsl = 1.3, close to the value at which the pair ampli- 
tude disappears in the elongated traps with A = 5 and 50. 

To show that the non-increasing behavior of Pqq for in- 
creasing A is not an artifact of the BdG approximation or 
a finite k-p^as, we employ the real-space self-consistent T- 
matrix approximation 02211 (RSTA), by which strongly inter- 
acting fermions in an inhomogeneous potential can be treated 
with high accuracy. RSTA has been shown to reproduce the 
pseudo-gap phase in high-Tc superconductors 12211 and the 
superconductor- insulator transition in disordered diamond su- 
perconductors li23ll . 



At the unitarity limit {l/{k-pas) — > 0), in the normal phase 
we self-consistently solve the following set of equations: 

Xsc(r,r') = r^Gt(r,r',w„)G4(r,r',-a;„) 



fir,r' 



Cir)d{r,r'), (5) 
[g~^S{r,r')+xsc{r,r')]~\ (6) 



S,(r, r', w„) = TG^ir', r, -uJn)f{r, r'), (7) 
G^ir,r',Lu„) = [[Gl}(r,r',w„)]"'-I],(r,r',w„)]"\?) 

where xsc is the pairing susceptibility, Ga (G°) the (non- 
interacting) Green's function, Y^a- the self-energy, w„ = 
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{2n + l)7r/r the Matsubara frequencies, and C(r) the space- 
dependent regularization factor, which is obtained as 



C(r; -^Max) = TT"^ / arctan i^jdk. 

(9) 

We discretize the system and use the rotational symmetry 
of the system to use a Fourier component expression in the 
relative azimuthal angle between two spatial lattice points. 
We need about 15 (positive) Matsubara frequencies, E'Max ^ 
30huj and 30 — 60 Fourier components for convergence at 
/i^ = lOSoJ > m and A = 1 or 12. 

If the phase transition from a normal gas to a superfluid is 
due to the divergence of the T-matrix, the maximum eigen- 
value Aiviax of xsc reaches zero from below at the transition 
point Tc- While in our trapped, finite-size system a first-order 
transition may occur, and then AMax is not necessarily zero, 
we believe that the transition should happen at similar values 
of Aiviax close to zero regardless of the trap aspect ratio A if 
the total number N and temperature T are similar Therefore, 
we compare AMax as a function of fi^ for a fixed fj,^ and T. 

As shown in Fig. |4] for A — 12, the value of Ai^ax is close 
to, but does not exceed, that for A = 1. This comparison is 
for aspect ratios smaller than those of MIT and Rice; however, 
because the effects of the trap shape are enhanced for smaller 
N, this result indicates that in the equilibrium, for ^ 10** 
atoms, the transition temperatures for a given P for A = 50 
should not exceed that for A = 5, and strengthens our conclu- 
sion that Pec is not enhanced as A is increased. 

To summarize, we have studied superfluidity of population- 
imbalanced fermions trapped in an axisymmetric harmonic 
trap by means of the Bogoliubov-de Gennes method. Our 
numerical results reproduce the major features of the experi- 
ments conducted at MIT, but does not reproduce those at Rice, 
as to the value of Pqc and as to the degree of LDA break- 
down. Recently Nascimbene and coworkers at ENS Paris 
II24II have trapped population-imbalanced fermions in elon- 
gated traps with various values of the aspect ratio and ob- 
served Pec = 0.76(3) and no deformation of density pro- 
files. Zwierlein and coworkers at MIT |25] observed very 
long spin diffusion time in two-species fermionic gases at uni- 
tarity by making two polarized gases collide in a quasi one- 
dimensional potential. The diffusion constant extracted from 
their experiment suggests that the timescale of the equilibra- 
tion is as long as one second for the configuration of the Rice 
experiment, which is much longer than the waiting time after 
the potential ramp in that experiment. We speculate that in 
the Rice experiment for P > 0.8, a non-equilibrium conden- 
sate state, possibly from the mechanism discussed in [8J, was 
observed in the course of slow relaxation and cooling process 
into a normal state. 
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